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Abstract. 

A theoretical model for the effect of water hydrogen bonding on the thermodynamics of hydrophobic 
hydration is proposed as a combination of the classical density functional theory with the recently developed 
probabilistic approach to water hydrogen bonding in the vicinity of a hydrophobic surface. The former allows 
one to determine the distribution of water molecules in the vicinity of a macroscopic hydrophobic particle and 
calculate the thermodynamic quantities of hydrophobic hydration as well as their temperature dependence, 
whereas the latter allows one to implement the effect of the hydrogen bonding ability of water molecules on 
their interaction with the hydrophobic surface into the DFT formalism. This effect arises because the number 
and energy of hydrogen bonds that a water molecule forms near a hydrophobic surface differ from their bulk 
values. Such an alteration gives rise to a hydrogen bond contribution to the external potential field whereto 
a water molecule is subjected in that vicinity. This contribution is shown to play a dominant role in the 
interaction of a water molecule with the surface. Our approach predicts that the free energy of hydration 
of a planar hydrophobic surface in a model liquid water decreases with increasing temperature in the range 
from 293 K to 333 K. This result is indirectly supported by the counter- intuitive experimental observation 
that under some conditions the hydration of a molecular hydrophobe is entropically favorable as well as by 
the molecular dynamics simulations predicting positive hydration entropy for sufficiently large (nanoscale) 
hydrophobes. 
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1 Introduction 

Particles exhibiting resistance to being either wetted by liquid water (case of meso/macroscopic 
ones) or dissolved therein (case of molecules) are generally referred to as hydrophobic solutes. The 
transfer of a hydrophobic solute into liquid water is accompanied by an increase in the free energy of 
the system, which results from the structural modifications of liquid water around the solute. This 
phenomenon is referred to as hydrophobic hydration. When two solute particles are close enough to 
each other, the total volume of water thus affected by both particles is smaller than when they are far 
apart. This gives rise to an effective, solvent-mediated (often referred to as hydrophobic) attraction 
between particles. 

It is believed that hydrogen bonding between neighboring water molecules 1-3 constitutes a key 
element of hydrophobic effects 4-7 (hydration and attraction) which, in turn, play a crucial role in 
many physical, chemical, and biological phenomena. 7-11 The development of predictive models ca- 
pable of estimating the temperature and pressure dependence of hydrophobic phenomena is therefore 
quite important. 

In attempts to understand hydrophobic effects at a fundamental level and to develop a generally 
satisfactory theory of hydrophobicity, 12,13 various mechanisms have been suggested 4-7 ' 14-21 (mostly 
involving the hydrogen bonding ability of water molecules). Despite many remaining controversies, 
the dependence of hydrophobic phenomena on the length scales of hydrophobic particles involved 
appears to be out of contention. 22-25 

Small hydrophobic molecules (whereof the linear size is comparable to that of a water molecule) 
can fit into the water hydrogen-bond network without destroying any bonds. 15 One one hand, this 
results in a negligible enthalpy of hydration; on the other hand, the presence of the solute is thought 
to constrain some degrees of freedom of the neighboring water molecules which should give rise 
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to a negative entropy of hydration proportional to the solute excluded volume. Consequently, the 
hydration free energy is positive and increases with temperature and solute excluded volume. 

Since the hydration of small hydrophobic molecules is entropically "driven" , so is their solvent- 
mediatated interaction. 12 ' 13 At small enough separations, the two hydrophobic molecules affect fewer 
solvent molecules than when they are far apart. Therefore, bringing two hydrophobic molecules 
sufficiently close enough to each other should result in a positive change of the entropy and should 
thereby lower the free energy of the solution (small enthalpy changes are neglected). 

Although attractive owing to its simplicity, such a mechanism of small scale hydrophobicity 
appears to be somewhat inaccurate. 12 ' 13 Simulations 26 ' 27 and theory 16 showed that two inert gas 
molecules (such as argon) would not be driven together to form a dimer; a solvent-separated pair 
would be a more likely state than a contact pair. This leads to a surprising suggestion that the 
hydration of small hydrophobic molecules is actually entropically favorable (the entropy of the system 
increases) which contradicts the conventional wisdom. 

The hydration of large hydrophobic particles is believed to occur via a different mechanism. 12,13,27 ~ : 
When inserted into liquid water, such a particle breaks some hydrogen bonds in its immediate vicinity. 
This results in a large positive enthalpy of hydration and hence in a free energy change proportional to 
the solute surface area (as opposed to being proportional to the solute volume for small hydrophobes). 

Thus, in contrast to entropically driven small-scale hydrophobicity, the hydration of large hy- 
drophobic particles is expected to be enthalpically driven and so is their hydrophobic interaction. 
Fewer water hydrogen bonds have to be broken when two large hydrophobes are "in contact" than 
when they are far from each other, so there is a negative enthalpy change when such particles aproach 
each other from larger separations. The free energy change (dominated by the enthalpy change) will 
be hence negative and will constitute a thermodynamic driving force for their attraction. 

Among theoretical means for studying hydrophobic phenomena (as well as many others including 
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phase transitions and phase equilibria), the methods of density functional theory 30,31 (DFT) have 
been particularly efficient. The DFT formalism has been widely used for studying the density profiles 
and thermodynamic properties of fluids near rigid surfaces of various sizes, shapes, and nature. 32 ' 33 

In DFT the interaction of fluid molecules with a foreign surface is usually treated in the mean- 
field approximation: every fluid molecule is considered to be subjected to an external potential 
that arises due to its pairwise interactions with the molecules of the impenetrable substrate. This 
external potential gives rise to a specific contribution to the free energy functional. The minimization 
of the free energy functional with respect to the number density of fluid molecules (as a function 
of the spatial coordinate r) provides their equilibrium spatial distribution. However, the effect of 
the impenetrable surface on the ability of fluid (water) molecules to form hydrogen bonds near the 
surface had been ignored so far in the conventional DFT formalism. 

In the present work (that can be considered as a sequel of our recently published paper 34 ), we at- 
tempt to fill in this gap and clarify some issues concerning the hydration of large hydrophobic particles 
by combining the density functional theory with the recently developed probabilistic appraoch 35 ' 36 
to hydrogen bonding between water molecules in the vicinity of a foreign surface. This approach 
provides an analytic expression for the average number of hydrogen bonds that a water molecule can 
form as a function of its distance to the surface. Knowing this expression, one can implement the 
effect of the hydrogen bonding of water molecules on their interaction with the hydrophobic surface 
into DFT, which is then employed to determine the distribution of water molecules near a macro- 
scopic hydrophobic particle and to calculate the thermodynamic quantities of hydrophobic hydration 
and their temperature dependence. 
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2 The outline of a probabilistic approach to water— water hydrogen 
bonding near a hydrophobic surface 

Let us first briefly describe the probabilistic hydrogen bond (PHB) model 35,36 for the hydrogen 
bonding ability of water molecules. It considers a water molecule, whereof the location is determined 
by its center, to have four hydrogen-bonding (hb) arms (each capable of forming a single hydrogen 
bond) of rigid and symmetric (tetrahedral) configuration with the inter-arm angles a = 109.47°. 
Each hb-arm can adopt a continuum of orientations. For a water molecule to form a hydrogen bond 
with another molecule, it is necessary that the tip of any of its hb-arms coincide with the second 
molecule. The length of an hb-arm thus equals the length rj of a hydrogen bond. 

The hydrogen bond length 77 is assumed to be independent of whether the molecules are in the 
bulk or near a hydrophobic surface. The characteristic length rj of pairwise interactions between 
water and molecules constituting the substrate (flat and large enough to neglect edge effects, with 
its location determined by the loci of the centers of its outermost, surface molecules) plays a simple 
role in the hydrogen bond contribution to hydration or hydrophobic interaction. 35 ' 36 It (rj) only 
determines the reference point for measuring the distance between water molecule and substrate, so 
it will be set equal to 77. 

Consider a "boundary" water molecule (BWM) in the vicinity of a hydrophobic substrate S 
(immersed in liquid water) at a distance x therefrom. Such a molecule forms a smaller number of 
hydrogen bonds (hereafter referred to as "boundary hydrogen bonds") than in bulk water because 
the hydrophobic surface restricts the configurational space available to other water molecules that 
are necessary for a BWM to form hydrogen bonds. The actual number of hydrogen bonds, that a 
particular BWM can form, depends on both its location and its orientation. A probabilistic hydrogen 
bond model allows one to obtain an analytic expression for the average number of bonds that a BWM 
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can form as a function of its distance to the surface ( "average" with respect to all possible orientations 
of the water molecule). 35 ' 36 

Note that a boundary hydrogen bond (BHB), involving at least one boundary water molecule, 
may be slightly altered energetically compared to the bulk one. Such alteration is still a sub- 
ject of contention 37 as different authors suggest opposite effects, i.e., both enhancement 18 ' 38 and 
weakening 30 of the boundary hydrogen bonds. In the PHB approach, there is no restriction on the 
energy of a bulk (water- water) hydrogen bond, £5 < 0, so that the approach is valid independent of 
whether e& < e s or e& = e s , or e& > e s , where e s < is the energy of a BHB. 

2.1 The average number of BHB's per water molecule 

Let us choose a Cartesian coordinate system so that its x-axis is normal to the plate P located at 
x = 0. Denote the number of hydrogen bonds per bulk water molecule by and the average number 
of hydrogen bonds per BWM by n s . The latter is a function of distance x between the water molecule 
and the hydrophobic surface (assumed to be smooth on a molecular scale): n s = n s (x). If x is larger 
than 2r], the number of hydrogen bonds that the molecule can form is not affected by the presence 
of the surface. Therefore, n s (x) = for x > 2rj. On the other hand, the function n s (x) attains its 
minimum at the minimal distance between the water molecule and the plate, i.e., at x = rj, because 
at this distance the configurational space available for the neighboring water molecules (to form a 
bond with the selected one) is restricted (compared to the bulk water) by the plate the most. The 
layer of thickness i] from x = rj to x = 2r\ is referred to as the surface hydration layer (SHL). 
The function n s = n s (x) can be shown 35 ' 36 to have the form 

n s = kih + k 2 b\ + k 3 bl + k A b\. (1) 

where b\ is the probability that one of the hb-arms (of a bulk water molecule) can form a hydrogen 
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bond and k\ = k±(x), &2 = ^(x), k% = ks(x), and k^ = k^{x) are coefficient-functions that can be 
evaluated by using geometric considerations (with their dependence on the BWM orientations being 
averaged). The functions k±(x), ^(x), k$(x), and k±(x) are presented in Figure 2a of ref.36. They all 
become equal to 1 at x > 1 where eq.(l) reduces to its bulk analog, rib = b\ + b\ + b\ + bf (see ref.39). 
Since experimental data on rib (and even its temperature dependence) are readily available, the latter 
equation allows one to determine the probability &i as its positive solution satisfying the condition 
< b\ < 1. Thus, equation (1) provides an efficient pathway to n s as a function of x. Figure 
2b in ref.36 presents this function for a hydrophobic flat surface immersed in water at temperature 
T = 293.15 K, which corresponds to rib = 3.65 hence b\ = 0.963707. 

The above expression for n s takes into account the constraint that some orientations of the 
hb-arms of a BWM cannot lead to the formation of hydrogen bonds because of the proximity to 
the hydrophobic particle. The severity of this constraint depends on the distance of the BWM to 
the surface, hence the x-dependence of ki,k2,k^, and k±. It assumes that the intrinsic hydrogen- 
bonding ability of a BWM (i.e., the tetrahedral configuration of its hb-arms and their lengths and 
energies) are unaffected by its proximity to the hydrophobic surface so that the latter only restricts 
the configurational space available to other water molecules necessary for this BWM to form hydrogen 
bonds. 

2.2 Hydrogen bond contribution to the interaction of a water molecule with a 
hydrophobic plate 

Knowing the function n s (x), one can examine the effect of water hydrogen bonding on the hydration 
of hydrophobic (and even composite) particles as well as on their solvent-mediated interaction. For 
example, let us derive an expression for ?7 e h x b t = U^ t (x), water- water hydrogen bond contribution 
to the total external potential field C/ cxt = f7 ext (a;) whereto a water molecule is subjected in the 
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vicinity of a hydrophobic surface. The latter is needed for the application of DFT methods to the 
thermodynamics of hydrophobic phenomena. 

The hydrogen bond contribution U^ t is due to the deviation of n s from n& (see refs.35 and 36) 
and the (possible) deviation of e s from e&. It can be determined as 

U£ = U% (x) = e s (x)n s (x) - e b n b (r, < x < oo) . (2) 

The first term on the RHS of this equation represents the total energy of hydrogen bonds of a water 
molecule at a distance x from the surface, whereas the second term is the energy of its hydrogen 
bonds in bulk (i.e., at x — > oo. Note that the dependence of U^ t on x may be due not only to 
the function n s (x), but also to the x-dependence of the hydrogen bond energy in the vicinity of the 
hydrophobic surface, e s = e s (x). In the PHB model n s {x) = n b for x > 2r] hence it is reasonable to 
assume that e s (x) = e& for x > 2r\ as well. Thus, U^ t {x) is a very short-ranged function of x, such 
that U^{x) = for x > 2r]. 

3 The outline of the methods of density functional theory 

The effect of water-water hydrogen bonding on the density profile of (liquid) water molecules in the 
vicinity of a hydrophobic surface can be now examined by using DFT. 30-33 In this formalism, the 
grand thermodynamic potential Q, of a nonuniform single component fluid, subjected to an external 
potential c7 cxt , can be represented as a functional of the number density p(r) of fluid molecules 

n[p(r)] = j v dvf h {p{v))+ 1 -j J dvdv'p{v)p{v')4>^\v-v'\) + J drU^{v)p{v)-p J dvp{v), (3) 

where V is the volume of the system, p is the equilibrium chemical potential, and </> at (|r — r'|) is 
the attractive part of the interaction potential between two fluid molecules located at r and r'. In 
this expression, the contribution to the free energy due to the short range repulsive interactions (the 
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first term on the RHS of the equation) is modeled by the hard sphere free energy in a local density 
approxiamtion (LDA) , with fh at r being the Helmholtz free energy density of a hard sphere fluid of 
uniform density equal to p(r). The longer ranged attractive interactions are treated in a mean-field 
(van der Waals) approximation and represented by the second term on the RHS of eq.(3). (Note 
that the LDA neglects short-ranged correlations which leads to the absence of oscillations in the 
density profile of a fluid near a hard wall. Although more accurate, nonlocal approximations are also 
available, 31 ' 40 ' 41 we preferred the LDA to ensure the transparency of presentation and to put the 
emphasis on the idea of combining the DFT methods with the PHB model). 

In an open thermodynamic system under conditions of constant chemical potential p, volume 
V, and temperature T (grand canonical ensemble), the equilibrium density profile is obtained by 
minimizing the functional fi[jo(r)] with respect to p(r), i.e., by solving the Euler-Lagrange equation 
5Q/5p = 0, which takes the form 



where Ph{p) = df\{p)/dp is the chemical potential of the uniform reference (hard sphere) fluid of 
density p. The substitution of the equilibrium density profile p(r) into eq.(3) provides the grand 
thermodynamic potential of the fluid. 

As already mentioned, the term f7 ext on the RHS of eq.(4) (and the corresponding term on the 
RHS of eq.(3)) had been conventionally meant to represent the external potential exerted by all 
the molecules constituting the hydrophobic substrate on a fluid molecule. Various models for the 
external potential were designed to take into account pairwise interactions of a fluid molecule with the 
molecules of the substrate 32 ' 33 as well as the effect of the latter on the pairwise interactions between 
fluid molecules themselves. 42 The contribution of these (pairwise) effects into U cxt will be denoted by 
UeZ to distinguish it from the hydrogen bond contribution, C/ c h x b t . Thus, the overall external potential 




(4) 
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whereto a water molecule is subjected near a hydrophobic surface can be represented as 



C/ cxt (r) = [/-(r) + Or) 



(5) 



In a closed thermodynamic system with constant number of molecules N, volume V, and temper- 
ature T (canonical ensemble), the chemical potential p,, appearing in eqs.(3) and (4), is not known in 
advance. Instead, it plays the role of a Lagrange multiplier corresponding to the constraint of fixed 
number of molecules in the system: 



N 



/ dvp{r). 
Jv 



This equation can be used to determine the Lagrange multiplier p, i.e., chemical potential in the 
system, as follows (see, e.g., ref.43). Introducing the "configurational" part of the hard sphere 
chemical potential as ph(p(r)) = ph(p(r)) — fcsTlnp(r), one can rewrite equation (3) in the form 

P(r) = exp [/i - C/ext(r) - p(p(r')) - J dv' p(v')^ t {\v - r'|)] }. 

Integrating this equation over the volume of the system and using the constraint on N, one obtains 

p = k B T\nN -k B Tln (| dr exp [ - C/ ext (r) - Ji{p{v')) - J dr / p(r')^ Bt (|r - r'|)] , 

whereof the substitution into eq. (3) yields the Euler-Lagrange equation for the density profile in the 
canonical ensemble: 



p h (p(r)) = k B TlnN- k B Tln (^J drexp [ - C/ ext (r) - /2(p(r')) - j dr'p(r')^. 

- f/ cxt (r)- J drV(r>at(|r-r'|), 



r — r 



(6) 



The Helholtz free energy F of the canonical ensemble can then be obtained by substituting the 
solution of eq.(6) into the corresponding functional of p(r): 

F[p(r)} = J v drf h (p(r)) + 1 - J J dvdv' p(v)p(v')^ t (\v - r'|) + J dr f/ cxt (r)p(r), (7) 
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In the particular case of (fluid) water near a flat hydrophobic surface, one can use the planar 
symmetry of the system and choose the Cartesian coordinates so that the surface is located in the 
y — z plane at x = with the molecules of the fluid occupying the "half-space" x > 0. The eqilibrium 
density profile obtained from eqs.(4) or (6) is then a function of a single variable x, i.e., p(r) = p(x). 

4 Free energy of hydration and its temperature dependence 

In a canonical ensemble, the free energy AF of hydration of a hydrophobic particle can be determined 
as the difference 

AF=F-F , (8) 

where F and Fq are the Helmholtz free energies of the system (liquid water) with and without a 
hydrophobic particle therein, respectively. Likewise, the free energy of hydrophobic hydration in a 
grand canonical ensemble can be determined as 

An = Q - fi , (9) 

where fi and CIq are the values of the grand thermodynamic potential of the system (liquid water) 
with and without a hydrophobic particle therein, respectively. 

Knowing the free energy of hydrophobic hydration, one can find <3? s and <3? E , the entropic and 
energetic contributions to AF, as 

= —TAS = T(dAF/dT) VN , $ E = AE = (0(AF/T)/0(1/T)) V , N , (10) 

respectively, such that AF = <3? E + <3? S (in eq.(10) the subscripts of the partial derivatives indicate the 
thermodynamic variables held constant upon taking the derivatives). Clearly, for the decomposition 
of AF into energetic and entropic components it is necessary to know its temperature dependence. 
(Note that a similar decomposition can be carried out for AQ.) 
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In the combined PHB/DFT-based model presented above, the temperature dependence of AF 
contains a contribution from the temperature dependence of C/ e h x b t , hydrogen bond contribution to the 
overall external field exerted by the hydrophobic surface on water molecules in its vicinity. As clear 
from eq.(2), the dependence of U^ t on T is due to the temperature dependence of four quantities: 
n s , rib, £s> and e&. The functions e& = Q,(T) and rib = ni,(T) are either readily available or can be 
constructed on the basis of available data. On the other hand, n s is unambiguously related to rib, 
hence its dependence on T can be considered to be known as well. Finally, it is reasonable to assume 
that, whether in the bulk or in the surface hydration layer, the energy of a hydrogen bond depends 
on temperature in such a way that the ratio e s (r)/e&(T) is independent of T. One can thus consider 
U^ t to be a known function of not only x but also T, C/ e h x b t = U^ t (x, T). This allows one to numerically 
determine the temperature dependence of AF and to subsequently use interpolation procedure to 
find an accurate analytical fit thereof which then can be used in eq.(10). 

5 Numerical evaluations 

In order to illustrate the above model with numerical calculations, we have considered the hydration 
of a flat, macroscopic hydrophobic surface in a model liquid water. The pairwise interactions between 
two water molecules 1 and 2 were modeled with the Lennard-Jones (LJ) potential, 



where r\i is the distance between molecules 1 and 2, e ww is the energy parameter and d is the diameter 
of a model molecule. The parameter e ww was adjusted to be 3.79 x 1CT 14 , which differs from its values 
used in the computer simulations (Monte Carlo or molecular dynamics) of various water models, in 
which e ww ranges 36 from 5.31 x 1(T 15 erg (ST2 model) to 1.47 x 1(T 14 erg (SWM4-NDPmodel) to 
2.54 x 10~ 13 erg (SSD model). Such a modification was needed to ensure that the phase diagram of 
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model water more or less resembles that of real water. In the DFT formalism this can be achieved 
only by adjusting the single intermolecular potential describing water- water interactions , whereas in 
computer simulations water-water interactions are usually described by the combination of LJ and 
electrostatic potentials; hence the difference in the energy parameters of the respective LJ potentials. 
The parameter d of the LJ potential also has different values in different water models 36 in the 
range from 3.02 A (SSD model) to 3.18 A (SWM4-NDP model). On the other hand, the length 77 
of the hydrogen bond (i.e., the distance between the oxygen atoms of two hydrogen-bonded water 
molecules) is reported 36 to be about 2.98 A. Since cPs (of various water models) and 77 are so close 
to each other, we assumed for our model d ~ 77. 

To find the equilibrium density profile of (model) water molecules in the vicinity of the hydropho- 
bic surface, it is necessary to solve eq.(4) using, e.g., an iterative procedure. 33 Namely, the density 
profile pi(x)at the i th iteration is found from the previous one, pi-\{x) via 

Ph{Pi{x)) = p - dr'pi-i(a;)^at(|r - r'|) - f7 ext (z). (11) 
Jv 

A similar iterative procedure can be used to solve equation (6), 

p h ( Pi (x)) = k B TlnN - k B Tln [J drexp { [ - U cxt (x) ~ Hpi-A*')) ~ J *'pi-i(x')<M|r - r'|)] 

- U ext (x)- J dr'p^ix'^dr-r'l). (12) 

For the chemical potential Ph{p) = ^h(p, T) of a hard sphere fluid we have adopted the well-known 
Carnaham-Starling approximation 32 ' 33 ' 44 
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» h (p,T) = fc B r(ln(A 3 p) + ^ * + ; 

where £ = (7rd 3 /6)p and A = (/i 2 /27rmfcsT) 1 / 2 is the thermal de Broglie wavelength of a model 
molecule of mass m (with h and k B being Planck's and Boltzmann's constants, respectively). Since 
Ph is a single- valued (monotonically increasing) function of p, one can extract Pi(x) from the LHS of 
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eq.(ll) or (12) and continue iterations. The attractive part $> at of the pairwise water-water interac- 
tions was modeled by using a well-known perturbation scheme: 45 

-e ww (ri2<2 1 /6^) ; 
0ww, (ri2 > 2 1 / 6 r/). 

Figure 1 presents the typical behaviour of C/ ext and its components. The thick dashed curve shows 
c7 cxt , while the lower and upper thin solid curves are for U^™ and U^ t , respectively. The function 
U^ t (x) was modeled as suggested in refs.32 and 33, 

oo (x < 7]), 



(13) 

-e sw exp[-A sw (x - rj)] (x > rj), 



where the energy parameter e sw was taken to be equal to A; p e ww with k p = 2.1 and the inverse length 
parameter A sw was set equal to 1/rj (the ratio k p = e sw /e ww characterizes the degree of hydrophobicity 
of the surface as the energy of a water molecule attraction to the surface at the distance rj between 
them relative to — e ww ). In C/ e h x b t (x) (see eq.(2)), the x-dependence of e s was approximated by a linear 
function increasing from its minimum value of k h €b at x = rj (with k h = 1.1 corresponding to slightly 
enhanced hydrogen bonds 18 ' 38 for molecules closest to the hydrophobic surface) to its maximum 
(bulk) value e& for x > 2ij. This results in U^ t (x) = ebrib[(k h — (k h — l)(x — r])/rj)n s (x)/nb — 1] for 
rj < x < 2r] and ?7 e h x b t (x) = for x > 2rj. 

As clear from Figure 1, the hydrogen bonding contribution to the external potential has a repulsive 
character unlike the conventional pairwise contribution that has an attractive character (note also 
that the former dominates the latter in the most part of the range i] < x < 2rj). The repulsive 
character of U^ t arises because that the total energy of hydrogen bonds per molecule near the 
hydrophobic surface is smaller (in absolute value) than in the bulk, which in turn is due to n s (x) < rib 
and e s (x) w e& for any x. 

At a given temperature, the evaluation of the free energy of hydration is simpler in a grand 
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canonical esemble, i.e., by solving eq.(ll), substituting the resulting equilibrium density profile in 
eq.(3), and then calculating AS1 according to eq.(9). This procedure was applied to the hydration 
of an infinitely large flat hydrophobic surface in the model liduid water at temperature T = 293.15 
K and chemical potential p = —11.5989 k B T corresponding to the vapor-liquid equilibrium of the 
model fluid. The liquid state of the bulk water was ensured by imposing the appropriate boundary 
condition p(x) — > p\ as x — > oo onto eq.(ll), with p\ being the bulk liquid density. The densities 
p v and pi of coexisting vapor and liquid, respectively, were found by solving a pair of equations 
expressing the conditions of phase equilibrium at a given T, 

P(P,T)\ p=Pv = KP,T)\ p=pi , p(p,T)\ p=Pv = P (p,T)\ p=pi , 

where 34 p(p, T) = pt(p, T)—ap and p(p, T) = pf l (p,T) — ^ap 2 and positive constant a = — J dr </> at (r). 
The pressure of a uniform hard sphere fluid ph = p^{p,T) is related to the corresponding chemical 
potential ph via dph/dp = pdph/dp. In the Carnahan-Starling approximation 34 ' 35 ' 41 

1 + £ + £2 _ £3 

Ph (p,T) = pk B T (1 _ g)3 • 

To clarify the effect of the two different contributions to U ext (x) on the water density distribution 
near the hydrophobe, the density profiles were obtained 34 by solving eq.(ll) with the overall external 
potential a) including both pairwise and hydrogen bond contributions (i.e., U ext (x) = U^ t (x) + 
U^ t (x)) and b) including only the pairwise component (i.e., with C/ cxt (x) = U^ t (x)). In both cases 
the ratio k p was taken to be 2.1. It was clearly demonstrated 36 that the hydrogen bond contribution 
to the external potential plays a crucial role in the formation of a thin, "strong depletion" layer (of 
density much lower than liquid and of thickness of a molecular diameter, in agreement with previous 
suggestions 12 ' 13 ) between liquid water and hydrophobic surface even for weakly hydrophobic surfaces 
(with high k p ). It was also shown that, as expected, even for a relatively strong hydrophobic surface 
(with low k p ) the conventional contribution to the external potential (due to pairwise interactions 
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between a water molecule and those of the substrate) cannot cause the formation of a vapor-like layer 
near the surface, although it does lead to a weak decrease in the vicinal fluid denstity compared to 
the bulk one. 

Figure 2 presents the grand canonical free energy of hydrophobic hydration AQ, expressed in 
units of ksT per rj 2 , as a function of the energetic alteration ratio k h of hydrogen bonds (in the SHL 
compared to the bulk) , at a constant ratio k p (characterizing the degree of the hydrophobicity of the 
surface). Due to the model character of the hydrophobic surface, we considered several values of k p 
(each curve in Figure 2a corresponds to a constant k p = 2.0,2.1,2.2 and 2.3 from bottom to top). 
As clear, the hydrophobic hydration is hardly sensitive to the hydrogen bond energy alteration ratio, 
k b , but quite sensitive to the degree of hydrophobicity of the surface, k p . The latter effect is also 
demonstrated in Figure 2b, where the dependence of AO on k p is shown for the case where in eqs.(3) 
and (11) only the pairwise component was included in the external potential, i.e., U ext (x) = U^ t (x). 

The temperature effect on the strength of solvent-mediated part of hydrophobic attraction is 
demonstrated by Figure 3 that presents the Helmholtz free energy of hydrophobic hydration, AF, 
and its energetic and entropic components, $ E and <£ s , as functions of T. The solid curve represents 
AF itself, while the long-dashed and short-dashed curves are for <£ E and <3? s , respectively. All the 
potentials are in units of ksT per rj 2 . The results in Figure 3a correspond to the overall external 
potential in eqs.(7) and (12) that includes both the pairwise and hydrogen bond contributions (i.e., 
U ext (x) = U^ t (x) + U^ t (x)), whereas Figure 3b shows the results obtained by solving eq.(12) and (7) 
with only the pairwise component included in the external potential (i.e., with U ext (x) = U^ t (x)). 
Both Figures 3a and 3b are for a hydrophobic surface with k p = 2.1. As clear, the water hydrogen 
bonding plays an important role in the hydrophobic hydration by making the process thermodynam- 
ically more unfavorable (i.e., significantly increasing the unfavorable free energy of hydration). 

As shown in Figure 3, the combined PHB/DFT model predicts the free energy of hydration of 
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a large hydrophobic particle to decrease with increasing temperature and suggests that the hydra- 
tion process is enthalpically unfavorable (i.e., the enthalpic contribution to the hydration free energy 
is positive), but entropically favorable (i.e., the entropic contribution to the hydration free energy 
is negative), with the latter effect being dominant. Currently, no experimental data on the ther- 
modynamics of hydration of large hydrophobes are available in literature. However, the enthalpic 
impediment to the hydration of a large hydrophobe is quite expected 12 ' 13,27-29 (being due to the 
breaking of vicinal hydrogen bonds). On the other hand, while the entropic enhancement of such 
hydration seems somewhat counter-intuitive, there are indirect experimental and simulational indica- 
tions of its physical soundness. For example, ref.44 reported experimental observation that dissolving 
an argon molecule in hot liquid water leads to an increase in entropy (i.e., the hydration of an argon 
molecule is entropically favorable), although the transfer of the same molecule into cold liquid wa- 
ter causes a decrease in entropy (i.e., its hydration is entropically unfavorable). These experiments 
are supported by a theoretical model (the two-dimensional Mercedes-Benz model with one fitting 
parameter) 46 as well as by the molecular dynamics simulations of SPC/E water model. 47 

Furthermore, studying the lengthscale dependence of hydrophobic hydration (under various ther- 
modynamic conditions) by means of MD simulations of SPC/E water model, it was demonstrated 23 
that the hydration thermodynamics changes its character from "entropy dominated" to "enthalpy 
dominated" near the crossover region as the length scale of a hydrophobe increases. At T = 300 K 
and pressure —1000 atm, the crossover region was found to be around R = 3 nm (R being the radius 
of a spherical hydrophobe). As reported, 23 the hydration is predominantly entropic (<fr s > $ E > 0), 
for solutes of radii smaller than 3 nm, whereas it is predominantly enthalpic (<5 E > <1> S > 0) for 
solutes of radii larger than 3 nm. As the radius of the solute increases, the enthalpic contribution 
<3? E increases whereas the entropic contribution $ s decreases and is expected to become negative "for 
sufficiently large solutes" (see ref.23). Qualitatively, this can be interpreted as the result of breaking 
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the tetrahedrally ordered structure of water-water hydrogen bond network by the foreign hydropho- 
bic particle; breaking the ordered structure of the hydrogen bond network is equivalent to increasing 
the disorder in the system which leads to an increase in its entropy. 

6 Concluding remarks 

In order to clarify some aspects of the effect of water-water hydrogen bonding on the thermodynamics 
of hydrophobic hydration, we have proposed a combination of our previously developed probabilistic 
approach to water-water hydrogen bonding with the classical density functional theory. The latter 
allows one to accurately determine the distribution of water molecules in the vicinity of a hydropho- 
bic particle and calculate the thermodynamic quantities of hydrophobic hydration as well as their 
temperature dependence. The former allows one to implement the effect of the hydrogen bonding 
ability of water molecules on their interaction with the hydrophobic surface into the DFT formalism. 

The hydrogen bond network of water molecules affects their interaction with the hydrophobic 
surface because the number and energy of hydrogen bonds that a water molecule forms in the the 
surface differ from their bulk values. Such an alteration gives rise to a short-range hydrogen bond 
contribution to the external potential field whereto a water molecule is subjected in that vicinity. 
This contribution is a dominant component of the interactions of a water molecule with the surface at 
distances between one and two hydrogen bond length. As we previously showed, 34 it plays a crucial 
role in the formation of a thin depletion layer (of thickness of a molecular diameter and of very low 
density, in agreement with previous suggestions 12,13 ) between liquid water and hydrophobic surface. 

The combined PHB/DFT approach to hydrophobic hydration predicts that the free energy of 
hydration of a model hydrophobic surface in a model liquid water decreases with increasing tem- 
perature in the range from 293 K to 333 K. It also corroborates the counter-intuitive experimental, 
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simulational, and theoretical observation that under some thermodynamic conditions the hydropho- 
bic hydration may be entropically favorable. As a possible explanation, one can conjecture that the 
destruction of the tetrahedrally ordered structure of water hydrogen bond network by a hydrophobic 
particle results in an increased disorder in the system which leads to an increase in its entropy. 

Acknowledgement - The author thanks Dr. E. Ruckenstein and Dr. G. Berim for many helpful discus- 
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Captions 

to Figures 1 to 3 of the manuscript "The effect of water-water hydrogen bonding on 

HYDROPHOBIC HYDRATION AND ITS TEMPERATURE DEPENDENCE" by Y. S. Djikaev and E. 

Ruckenstein. 

Figure 1. The typical behaviour of the overall external potential U axt (exerted by a hydrophobic 
surface on water molecules in its vicinity), shown as a thick dash-dotted curve, and its components, 
[/ e p x * (lower thin solid curve) and ?7 B h x b t (upper thin solid curve). 

Figure 2. The grand canonical free energy of hydrophobic hydration AQ, expressed in units of ksT 
per rj 2 : a) as a function of the bond energy alteration ratio k h at a constant k p . Different curves 
correspond to different degrees of hydrophobicity {k p = 2.0, 2.1, 2.2, and 2.3 from bottom to top); 
b) as a function of k p for the case where U ext (x) = U^(x). 

Figure 3. The Helmholtz free energy of hydrophobic hydration, AF, and its energetic and entropic 
components, <3? E and $ s , as functions of T for a hydrophobic surface with k p = 2.1 and a) U cxt (x) = 
U eZ( x )) and b ) U ext (x) = U»™(x)+U^(x). The solid curves represent AF itself, while the long-dashed 
and short-dashed curves are for <3> E and $ s , respectively. 
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